In [1]:
cd C:\Users\alpha\Documents\work\burned pixels\ecoregions_countries_climates_continents

C:\Users\alpha\Documents\work\burned pixels\ecoregions_countries_climates_continents
In [2]:
import geopandas as gpd
import matplotlib.pyplot as plt
import pickle
import pandas as pd
import numpy as np
import calendar

import matplotlib.cm as cm
import matplotlib.colors as cls
from mpl_toolkits.axes_grid1.inset_locator import inset_axes
In [3]:
#latex
import matplotlib
matplotlib.rcParams['mathtext.fontset'] = 'custom'
matplotlib.rcParams['mathtext.rm'] = 'Times New Roman'
matplotlib.rcParams['mathtext.it'] = 'Times New Roman:italic'
matplotlib.rcParams['mathtext.bf'] = 'Times New Roman:bold'
In [4]:
from scipy.stats import boxcox, probplot
In [5]:
bios = gpd.read_file('biomes_hemispheres')
In [6]:
bios.columns
Out[6]:
Index(['OBJECTID', 'BIOME_NUM', 'is_north', 'SUM_2001_J', 'SUM_2001_F',
       'SUM_2001_M', 'SUM_2001_A', 'SUM_2001_1', 'SUM_2001_2', 'SUM_2001_3',
       ...
       'SUM_2017_3', 'SUM_2017_4', 'SUM_2017_S', 'SUM_2017_O', 'SUM_2017_N',
       'SUM_2017_D', 'SUM_Area', 'Shape_Leng', 'Shape_Area', 'geometry'],
      dtype='object', length=211)
In [7]:
area = bios['SUM_Area']*1e-6
In [10]:
data_columns = bios.columns[3:3+17*12]
In [11]:
inter_cols = {}
for yr in range(17):
    inter_cols[yr] = data_columns[12*yr:12*yr +12]
In [12]:
intra_cols = {}
for mo in range(12):
    intra_cols[mo] = [data_columns[mo + yr *12  ] for yr in range(17)]
In [13]:
for yr in range(17):
    inter = bios[  inter_cols[yr]].sum(axis=1)/area
    label = 'inter_' + str(yr +2001)
    bios= bios.assign(label=inter)
    bios = bios.rename(columns={'label': label})
In [14]:
inter_columns = bios.columns[-17:]
In [15]:
inter_max = bios[inter_columns].max().max()
In [16]:
for mo in range(12):
    intra = bios[ intra_cols[mo]].sum(axis=1)/(17*area)
    label = 'intra_{:02}'.format(mo)
    bios = bios.assign(label=intra)
    bios = bios.rename(columns={'label': label})
In [17]:
intra_columns = bios.columns[-12:]
In [18]:
intra_max = bios[intra_columns].max().max()
In [62]:
def draw_map(df,column,title=None,vmin=0,vmax=1,name='name',dpi=100,cmap='OrRd'):
    fig, ax = plt.subplots(figsize=(12, 8))

    ax.set_aspect('equal')
    ax = df.plot(ax =ax, column=column, cmap=cmap,vmin=vmin, vmax=vmax)
    ax.set_axis_off()
    
    norm = cls.Normalize(vmin, vmax)
    cmmapable = cm.ScalarMappable(norm, cmap)
    
    cmmapable._A = []
    cbaxes = inset_axes(ax, width="80%", height="1%", loc=3)
    cb = plt.colorbar(cmmapable, cax=cbaxes,orientation="horizontal") 

    cb.set_label(title, fontsize=20, family='Times New Roman')
    cb.ax.set_yticklabels(cb.ax.get_yticklabels(), fontsize=15, family='Times New Roman')

    plt.show()
    fig.savefig(name, dpi=dpi, bbox_inches='tight')
In [ ]:
''' backup
def draw_map(df,column,title=None,vmin=0,vmax=1,name='name',dpi=100):
    fig, ax = plt.subplots(figsize=(12, 8))

    ax.set_aspect('equal')
    ax = df.plot(ax =ax, column=column, cmap='OrRd',vmin=vmin, vmax=vmax)
    ax.set_axis_off()
    
    norm = cls.Normalize(vmin, vmax)
    cmmapable = cm.ScalarMappable(norm, 'OrRd')
    
    cmmapable._A = []
    cbaxes = inset_axes(ax, width="80%", height="1%", loc=3)
    cb = plt.colorbar(cmmapable, cax=cbaxes,orientation="horizontal") 

    cb.set_label(title, fontsize=20, family='Times New Roman')
    cb.ax.set_yticklabels(cb.ax.get_yticklabels(), fontsize=15, family='Times New Roman')

    plt.show()
    fig.savefig(name, dpi=dpi, bbox_inches='tight')
'''

Inter-anual

In [20]:
lambdas = []
for yr in range(17):
    column= 'inter_{}'.format(2001 + yr)
    test_column, lbd = boxcox(bios[column] + 1e-25)
    lambdas.append(lbd)
lmbda=  sum(lambdas)/len(lambdas)
In [21]:
lambdas
Out[21]:
[0.11549818404472534,
 0.13247022652837614,
 0.14040360802627982,
 0.1317583680617727,
 0.13043098141653348,
 0.13325660542853673,
 0.12843417474371077,
 0.12925781778553824,
 0.1306853957846011,
 0.11861652800316695,
 0.12884257581501873,
 0.13166418002449554,
 0.13013977100347757,
 0.11986988692310972,
 0.12033228362428698,
 0.1314218302256409,
 0.1224914409660524]
In [22]:
probplot(bios.inter_2017, plot=plt.gca());
In [23]:
probplot(test_column, plot=plt.gca());
In [24]:
maxs = []
mins = []
for yr in range(17):
    column= 'inter_{}'.format(2001 + yr)
    test_column = boxcox(bios[column] + 1e-25,lmbda=lmbda)
    maxs.append(test_column.max())
    mins.append(test_column.min())

box_min = min(mins)
box_delta =max(maxs) - box_min
In [25]:
for yr in range(17):
    column = 'inter_{}'.format(yr +2001)
    test_column = (boxcox(bios[column] + 1e-25,lmbda=lmbda) - box_min)/box_delta
    bios = bios.assign(test_column =test_column)
    yr_name = str(yr +2001)
    #"Burnt area pixels (#/km2) by biome in 2001
    title = r'Burned area index (#/$\mathrm{km}^{2}$) by biome N-S parts in ' +yr_name
    #title = r'Burned Area Pixels/$\mathrm{Km}^{2}$ by Climate 2001-2017 {} '.format(calendar.month_name[mo + 1])
    #title =(r'ABC123 vs $\mathrm{ABC123}^{123}$')
    name = 'inter_biomes_nsparts_index_{}'.format(yr +2001)
    draw_map(bios,'test_column',title=title,name=name)

C:\ProgramData\Anaconda3\envs\geopandas\lib\site-packages\matplotlib\font_manager.py:1339: UserWarning: findfont: Font family ['cursive'] not found. Falling back to DejaVu Sans
  (prop.get_family(), self.defaultFamily[fontext]))

Intra-anual

In [26]:
lambdas = []
for mo in range(12):
    column= 'intra_{:02}'.format(mo)
    test_column, lbd = boxcox(bios[column] + 1e-25)
    lambdas.append(lbd)
lmbda=  sum(lambdas)/len(lambdas)
In [27]:
lambdas
Out[27]:
[0.09657411870392897,
 0.09397366386268413,
 0.12691276070279361,
 0.11853975794400545,
 0.13844976415288773,
 0.1244002609968398,
 0.11776589776775662,
 0.11683053353063909,
 0.11450503528246088,
 0.10609081953392038,
 0.10146759197041494,
 0.09515963130263318]
In [28]:
probplot(bios.intra_11, plot=plt.gca());
In [29]:
probplot(test_column, plot=plt.gca());
In [30]:
maxs = []
mins = []
for mo in range(12):
    column= column= 'intra_{:02}'.format(mo)
    test_column = boxcox(bios[column] + 1e-25,lmbda=lmbda)
    maxs.append(test_column.max())
    mins.append(test_column.min())

box_min = min(mins)
box_delta =max(maxs) - box_min
In [31]:
for mo in range(12):
    column = 'intra_{:02}'.format(mo)
    test_column = (boxcox(bios[column] + 1e-25,lmbda=lmbda) - box_min)/box_delta
    bios= bios.assign(test_column =test_column)
    mo_name = calendar.month_name[mo + 1]
    #"Burnt area pixels (#/km2) by biome in Janury (2001-2017)"
    title = r'Burned area index (#/$\mathrm{km}^{2}$) by biome N-S parts in ' +mo_name + ' (2001-2017)'
    #title = r' Burned Area Pixels/$\mathrm{Km}^{2}$ by climate 2001-2017 ' + mo_name
    #title = r'Burned Area Pixels/$\mathrm{Km}^{2}$ by Climate 2001-2017 {} '.format(calendar.month_name[mo + 1])
    #title =(r'ABC123 vs $\mathrm{ABC123}^{123}$')
    name = 'intra_biomes_nsparts_index_{:02}'.format(mo)
    draw_map(bios,'test_column',title=title,name=name)

Total

In [32]:
total = bios[data_columns].sum(axis=1)/17

total = total/area

bios= bios.assign(total=total)
In [33]:
column = 'total'
test_column, lbd = boxcox(bios[column] + 1e-25)
test_column = (test_column - test_column.min())/(test_column.max() - test_column.min())
bios = bios.assign(test_column =test_column)
title = r'Burned area index (#/$\mathrm{km}^{2}$) by biome N-S parts in (2001-2017)' 
#title = r'Burned Area Pixels/$\mathrm{Km}^{2}$ by Climate 2001-2017 {} '.format(calendar.month_name[mo + 1])
#title =(r'ABC123 vs $\mathrm{ABC123}^{123}$')
name = 'total_biomes_nsparts_index'
draw_map(bios,'test_column',title=title,name=name)

Anomalias

Inter-anual

In [38]:
inter_medias = bios[inter_columns].mean(axis=1)
In [41]:
matriz = np.array([[4,6,8],[6,8,10]])
In [43]:
matriz - np.array([[0],[2]])
Out[43]:
array([[4, 6, 8],
       [4, 6, 8]])
In [53]:
matriz.shape
Out[53]:
(2, 3)
In [54]:
np.array([[0],[2]]).shape
Out[54]:
(2, 1)
In [58]:
anomalias = np.array(bios[inter_columns]) -  inter_medias.values.reshape(27,1)
In [45]:
bios[inter_columns]
Out[45]:
inter_2001 inter_2002 inter_2003 inter_2004 inter_2005 inter_2006 inter_2007 inter_2008 inter_2009 inter_2010 inter_2011 inter_2012 inter_2013 inter_2014 inter_2015 inter_2016 inter_2017
0 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000e+00 0.000000 0.000000 0.000000 0.000000
1 0.000000 0.000005 0.000010 0.000003 0.000004 0.000016 0.000003 0.000003 0.000031 0.000000 0.000003 0.000007 5.462977e-07 0.000000 0.000000 0.000001 0.000000
2 0.047521 0.076068 0.076465 0.079402 0.093502 0.071258 0.096334 0.061346 0.064016 0.093142 0.058514 0.064990 5.138119e-02 0.058642 0.074753 0.058978 0.060809
3 0.034023 0.039508 0.073047 0.083512 0.072145 0.058385 0.087838 0.056454 0.074407 0.082565 0.052387 0.067192 6.126097e-02 0.062217 0.070823 0.070867 0.055229
4 0.072870 0.087703 0.084904 0.122496 0.116345 0.082425 0.149634 0.081923 0.060723 0.150865 0.069885 0.100991 6.301827e-02 0.061261 0.105633 0.098319 0.076096
5 0.086836 0.094555 0.125449 0.154393 0.106404 0.105181 0.119554 0.112571 0.116948 0.098936 0.123857 0.118240 1.001687e-01 0.093120 0.113729 0.091015 0.121322
6 0.035056 0.042468 0.069894 0.027477 0.089520 0.050962 0.043930 0.061433 0.054920 0.030546 0.131657 0.068319 5.678596e-02 0.026642 0.020432 0.048741 0.069855
7 0.029098 0.048656 0.054778 0.024103 0.016866 0.064001 0.017241 0.008843 0.031539 0.006307 0.041213 0.026789 3.237772e-02 0.024398 0.018662 0.024845 0.033752
8 0.021624 0.028188 0.026721 0.027089 0.031283 0.038570 0.030720 0.050702 0.041731 0.035059 0.029887 0.033241 2.128196e-02 0.035199 0.025436 0.017554 0.021997
9 0.005998 0.026541 0.025904 0.005555 0.008045 0.012225 0.018509 0.027902 0.009625 0.009517 0.016561 0.026370 1.761131e-02 0.015946 0.020750 0.010911 0.027393
10 0.010041 0.024139 0.037167 0.014423 0.014518 0.015511 0.009663 0.019393 0.011776 0.013669 0.018976 0.029045 2.327012e-02 0.027830 0.021612 0.016264 0.015823
11 0.748677 0.696804 0.721054 0.822319 0.739055 0.723421 0.796940 0.707566 0.695598 0.754241 0.782219 0.828883 6.526381e-01 0.683454 0.684785 0.626272 0.668016
12 0.584928 0.592581 0.670168 0.620278 0.735867 0.589431 0.649654 0.613481 0.514674 0.528638 0.553869 0.509473 5.194914e-01 0.485283 0.468422 0.585141 0.473098
13 0.106718 0.046370 0.049552 0.019848 0.021067 0.015770 0.006652 0.015294 0.016340 0.006377 0.020399 0.016747 1.491790e-02 0.023587 0.016976 0.036117 0.050066
14 0.072921 0.116338 0.084369 0.074839 0.091262 0.097523 0.082140 0.112901 0.052278 0.076596 0.054199 0.060166 2.074878e-02 0.078761 0.074685 0.040418 0.073410
15 0.771776 0.964372 0.686142 0.812639 0.937373 0.757694 0.824067 0.797177 0.719120 0.847723 0.660443 0.784084 6.730329e-01 0.589442 0.725515 0.721452 0.643498
16 0.541217 1.063833 0.618616 0.795856 0.700235 0.703267 0.784351 0.793291 0.962006 0.442947 0.828018 0.729685 6.371421e-01 0.677844 0.946413 0.694031 0.635496
17 0.079299 0.090975 0.112876 0.069274 0.092837 0.080094 0.076181 0.078217 0.068832 0.090931 0.083040 0.069390 6.981568e-02 0.070018 0.063516 0.048201 0.066648
18 0.002930 0.006434 0.005321 0.003741 0.003398 0.004782 0.004421 0.003632 0.003080 0.003485 0.003976 0.003964 2.306099e-03 0.002411 0.003797 0.002171 0.002025
19 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000 0.000000e+00 0.000000 0.000000 0.000000 0.000000
20 0.001741 0.002741 0.009585 0.010173 0.003218 0.001153 0.002208 0.001646 0.001829 0.006873 0.001070 0.001654 1.318784e-03 0.000751 0.003388 0.004603 0.002067
21 0.041173 0.068758 0.043392 0.061789 0.058256 0.045308 0.024458 0.026639 0.046589 0.052355 0.026336 0.036716 1.943181e-02 0.051047 0.055112 0.036625 0.073661
22 0.019080 0.016105 0.034172 0.019516 0.019392 0.013835 0.034651 0.016257 0.024260 0.015301 0.014927 0.024276 1.792923e-02 0.014123 0.022806 0.021058 0.033771
23 0.372587 0.408419 0.045588 0.130526 0.044709 0.187594 0.159600 0.107300 0.042331 0.108560 0.506366 0.363277 7.190637e-02 0.119366 0.125376 0.033206 0.219797
24 0.009365 0.029262 0.015035 0.023997 0.025421 0.025155 0.017233 0.013051 0.013034 0.018543 0.014165 0.012830 9.764005e-03 0.013567 0.015172 0.011916 0.027799
25 0.010647 0.066837 0.040874 0.059064 0.042183 0.093326 0.029282 0.024870 0.046865 0.028810 0.043235 0.040294 4.656843e-02 0.054597 0.088643 0.039079 0.027083
26 0.031949 0.038983 0.086925 0.077045 0.070396 0.077495 0.069128 0.054948 0.076532 0.061753 0.064178 0.055275 6.154448e-02 0.075962 0.086954 0.070002 0.051241
In [61]:
anomalias[0,:]
Out[61]:
array([0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0.])
In [56]:
 inter_medias.values.reshape(27,1)
Out[56]:
array([[0.00000000e+00],
       [5.14162559e-06],
       [6.98307459e-02],
       [6.48152618e-02],
       [9.32406375e-02],
       [1.10722165e-01],
       [5.46259154e-02],
       [2.96158009e-02],
       [3.03695063e-02],
       [1.67860444e-02],
       [1.90069836e-02],
       [7.25408440e-01],
       [5.70263335e-01],
       [2.83998683e-02],
       [7.43267022e-02],
       [7.59738282e-01],
       [7.38485237e-01],
       [7.70673324e-02],
       [3.63958779e-03],
       [0.00000000e+00],
       [3.29509800e-03],
       [4.51556728e-02],
       [2.12623553e-02],
       [1.79206473e-01],
       [1.73711282e-02],
       [4.60151734e-02],
       [6.53124234e-02]])
In [64]:
anomalias.max()
Out[64]:
0.32715982717537395
In [65]:
for yr in range(17):
    test_column = anomalias[:,yr]
    bios = bios.assign(test_column =test_column)
    yr_name = str(yr +2001)
    #"Burnt area pixels (#/km2) by biome in 2001
    title = r'Anomalies burned pixels (#/$\mathrm{km}^{2}$) by biome N-S parts in ' +yr_name
    #title = r'Burned Area Pixels/$\mathrm{Km}^{2}$ by Climate 2001-2017 {} '.format(calendar.month_name[mo + 1])
    #title =(r'ABC123 vs $\mathrm{ABC123}^{123}$')
    name = 'inter_anomalies_biomes_nsparts_index_{}'.format(yr +2001)
    draw_map(bios,'test_column',title=title,name=name,cmap='seismic',vmin=anomalias.min(),vmax=anomalias.max())

Intra-anual

In [69]:
intra_medias = bios[intra_columns].mean(axis=1)
In [70]:
anomalias = np.array(bios[intra_columns]) -  intra_medias.values.reshape(27,1)
In [71]:
for mo in range(12):
    test_column = anomalias[:,mo]
    bios= bios.assign(test_column =test_column)
    mo_name = calendar.month_name[mo + 1]
    #"Burnt area pixels (#/km2) by biome in Janury (2001-2017)"
    title = r'Anomalies burned area pixels (#/$\mathrm{km}^{2}$) by biome N-S parts in ' +mo_name + ' (2001-2017)'
    #title = r' Burned Area Pixels/$\mathrm{Km}^{2}$ by climate 2001-2017 ' + mo_name
    #title = r'Burned Area Pixels/$\mathrm{Km}^{2}$ by Climate 2001-2017 {} '.format(calendar.month_name[mo + 1])
    #title =(r'ABC123 vs $\mathrm{ABC123}^{123}$')
    name = 'intra_anomalies_biomes_nsparts_index_{:02}'.format(mo)
    draw_map(bios,'test_column',title=title,name=name,cmap='seismic',vmin=anomalias.min(),vmax=anomalias.max())

Total

In [73]:
anomalia = bios.total -bios.total.mean()
In [74]:
anomalia
Out[74]:
0    -0.142369
1    -0.142364
2    -0.072538
3    -0.077554
4    -0.049128
5    -0.031647
6    -0.087743
7    -0.112753
8    -0.112000
9    -0.125583
10   -0.123362
11    0.583039
12    0.427894
13   -0.113969
14   -0.068042
15    0.617369
16    0.596116
17   -0.065302
18   -0.138729
19   -0.142369
20   -0.139074
21   -0.097213
22   -0.121107
23    0.036837
24   -0.124998
25   -0.096354
26   -0.077057
Name: total, dtype: float64
In [75]:
test_column = anomalia
bios= bios.assign(test_column =test_column)

title = r'Anomalies burned area pixels (#/$\mathrm{km}^{2}$) by biome N-S parts in (2001-2017)'

name = 'total_anomalies_biomes_nsparts'
draw_map(bios,'test_column',title=title,name=name,cmap='seismic',vmin=anomalia.min(),vmax=anomalia.max())